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ABSTRACT 

The solar wind emanates from the hot and tenuous solar corona. Earlier studies using 1.5 dimensional 
simulations show that Alfven waves generated in the photosphere play an important role in coronal heating 
through the process of non-linear mode conversion. In order to understand the physics of coronal heating and 
solar wind acceleration together, it is important to consider the regions from photosphere to interplanetary space 
as a single system. We performed 2.5 dimensional, self-consistent magnetohydrodynamic simulations, covering 
from the photosphere to the interplanetary space for the first time. We carefully set up the grid points with 
spherical coordinate to treat the Alfven waves in the atmosphere with huge density contrast, and successfully 
simulate the solar wind streaming out from the hot solar corona as a result of the surface convective motion. The 
footpoint motion excites Alfven waves along an open magnetic flux tube, and these waves traveling upwards in 
the non-uniform medium undergo wave reflection, nonlinear mode conversion from Alfven mode to slow mode, 
and turbulent cascade. These processes leads to the dissipation of Alfven waves and acceleration of the solar 
wind. It is found that the shock heating by the dissipation of the slow mode wave plays a fundamental role in 
the coronal heating process whereas the turbulent cascade and shock heating drive the solar wind. 

Subject headings: Sun: photosphere — Sun: chromosphere — Sun: corona — stars: mass-loss 



Introduction 



■ The coronal heating and solar wind acceleration arc 
"fundamental problems in solar physics. Although various 
physical mechanisms have been proposed for coronal heat- 
"ing, it remains unclear why the hot corona exists above the 
C ■ 'cool photosphere and the high-speed solar wind emanates 
.from there. The main difficulty arises due to the rapid 
-decrease of the density, amounting to more than 15 orders 
C^) of magnitude in the interplanetary space compared to the 
t-H ,photospheric value. The huge density contrast between 
t-H -the photosphere and interplanetary space has made the 
[problem difficult to understand the energy transfer from 
. — i -the Sun to the solar wind as a single system even within 
the magnetohydrodynamic (MHD) framework. 
Jjj Alfven wave is believed to be a primary candid ate that 
. . . -drives the solar wind (e.g. iMcIntosh et al.l f201 lh . Direct 
observations of propagating Alfven wave s have been re- 
ported after the launch of Hinode satellite (iDe Pontieu et all 
120071 : IRishizuka et al.l2008HOkamoto fc De Pontieul2011| ). 
Since Alfven waves are notoriously difficult to dissipate, 
various physical processes have been proposed. The dissi- 
pation of Alfven waves is the key behind the acceleration 
of solar wind and the essential problem is to dissipate the 
Alfven wave which eventually will transfer the energy to 
accelerate the solar wind. Recently, the mechanism of 
turbulent cascade has been proposed in which the down- 
ward wave that is generated due to the reflection of Alfven 
wave in the gravitationally stratified atmosphere interacts 
with the upward propagating Alfven wave and develops 
Alfvcnic turbulence. Once the Alfvenic turbulence is gen- 
erated, the energy cascade of the turbulence to smaller spa- 



tial scales finally heats t he corona and drives the solar wind 
(Mattha cus et al.l Il999f ). Various phenomenological ap- 
proaches have been considered to avoid the complexities of 
the turbulence (ICranmer et al.l l2007t IVerdini k Veliill2007t 
iBigot et al.l |2008|) . In numerical simulati ons an incom - 



pressible approximation is usually adopted dE inaudi et al 
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20111 ). although the mode conversion fr om Alfven to slow 
mode seems to play an important role dKudoh fc Shibata 
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lAntolin fc Shibata 
Matsumoto fc Shibatal I2010D . One of the major 
challenges in numerical simulation is to consider the huge 
density contrast between the solar photosphere and inter- 
planetary space. In this paper, we show results of two- 
dimensional MHD simulations, considering region from 
the solar photosphere and solar wind as a single system 
and include the details of wave reflection from the tran- 
sition region, nonlinear mode conversion as well as the 
turbulent cascade for the first time. 

2. Method 

We have performed two-dimensional MHD simulation 
with radiative cooling, thermal conduction, and grav- 
ity. Our numerical model includes a single flux tube ex- 
tended from a kilo G auss (kG) patch in the polar region 
(jTsuneta et al.l l2008j ) to the interplanetary space (~20 
Rs). Our basic equations are 



f + V.(pV)=0, 
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Fig. 1. — Results of MHD simulation of Alfven wave propagation from the solar photosphere to the interplanetary space, 
(a) Normalized density fluctuation, (b) region that is magnified 5 times. The red squared region in (a) is equivalent to (b). 
(c) Temperature distribution, (d) density distribution, the regions shown are magnified 100 times from (b). The red squared 
region in (b) is equivalent to (c) and (d). The white solid lines in each panel represent the magnetic field lines. Lengths are 
shown in units of Rs = 6.96 x 10 5 km. 

(r, 6 = 7r/2,</>) and dg = 0. Initial magnetic field is ex- 
trapolated using potential field approximation so that the 
open field lines are extended from a kilo Gauss strength 
to 10 Gauss at the height of 2,000 km. Initially, we 
set an isothermal (10 4 K) atmosphere. We input ve- 
locity disturbance in 9 direction at the footpoint of the 
flux tube in order to generate Alfven waves. We assume 
white noise power spectrum with (1 km s^ 1 ) 2 in total 
power, which is r ather a good approx i matio n to the ob- 
served spectrum (jMatsumoto fc Kital I2010D . Since the 
veloci ty disturbances are p re-determined, any feedback ef- 
fects (jGrappin et al.ll20 08) of the coronal disturbances on 
the photosphcric motions are ignored. Periodic bound- 
ary condition is posed in the <p directi on. Our numerical 
simu lation is based on HLLD scheme (jMivoshi fc Kusanol 
120051 ) that is robust and efficient among the various kind 
of approximate Ricmann solvers. The solcnoidal condi- 
tion ,V • B = 0, is sati sfied within a round-off error by 
using flux-CT method <[T6thl l2000h . The TVD-MUSCL 
scheme enables us to archive the second order accuracy 
in space, while the Runge-Kutta method gives the second 
order accuracy in time. We use min-mod limitcr in order 
to suppress the numerical oscillation around shocks, which 
is one of the standard technique in TVD-scheme. Our nu- 
merical domain extends from the photosphere to R = 20 
Rs radially. The spacial resolution is 25 km at the surface 
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where, p, V, B, and T are mass density, velocity, mag- 
netic field, and temperature, respectively Pt indicates 
total pressure, P g + B 2 /8tt, where P g is gas pressure. To- 
tal energy per unit volume is described as £ = pV 2 /2 + 
Pg/d — 1) + B 2 /8ir with specific heat ratio 7 = 5/3. g 
is gravitational acceleration, —GAlr/r 2 , where G and M 
are the gravitational constant and solar mass, respectively. 
We adopt anisotropi c thermal conduction tenso r k along 
magnetic field lines ([Yokovama fc Shibatall200l . Radia- 
tive cooling term R is as sumed to be a function o f local 
density and temperature (|Suzuki fc Inuts uka 2005|). Note 
that no other ad hoc source terms for heating are included 
in the energy equation. 

The vector form of our basic equations are appropri- 
ately transformed into spherical coordinate system with 
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and increasing with radius. The horizonal length is 3,000 
km at the photosphere with spacial resolution of ~ 100 
km. Total grid points in our simulation are 8198 in radial 
direction and 32 in horizontal direction. 

The turbulent heating rate is estimated by dimensional 
analysis since viscosity and resistivity are not included ex- 
plicitly in our basic equations. First, we derive the Fourier 
component, vg, of the Alfvenic disturbance (vg) as follows. 

vg{r,h)= J vg(r^)e- tk ^rdcp (5) 

Then, energy spectral density, E(r,k$), becomes 

E(r, M = [Mr, fc )| 2 + \v e (r, -fc )| 2 ] , (6) 

where A(j) indicates the angular system size in direction. 
By using and E(r, k^), we can estimate the energy ex- 
changing rate, e(r, fc^), for a certain wave number, fc^, with 
neighboring Fourier modes. Then e(r, k^) becomes 

e(r,fc )~p£(r,fc ) 3 / 2 fc;f , (7) 

where p denotes the mean density averaged over time and 
(j) direction. As a turbulent heating rate, we use e(r, k^) 
whose wave number is larger than a critical wave number 
that is determined by numerical resolution. We choose 
k ( f,rA(f>/2TT = 4 as the critical wave number that corre- 
sponds to the spatial resolution covering one wave length 
by 8 grid points in our simulation. 

3. Results &l Discussions 

The coupling between the magnetic field and the sur- 
face convection exc ites upward propagating Alfven wave 
(Stc iner et al.lll998l) . and it could be an efficient energy 
carrier in the solar atmosphere. Considering such scenario, 
we performed 2.5 D MHD simulations covering regions di- 
rectly from the photosphere to the interplanetary space. 
Once the Alfven wave is forced to excite, the numerical 
system attains quasi-steady state within 1,800 minutes. 
Due to the dissipation of the Alfven wave, the initially 
static and isothermal (10 4 K) atmospheres eventually de- 
velops a hot corona (10 6 K) and a high-speed (> 500 km 
s _1 ) solar wind (Figs. 1 and 2). The radial profiles of ve- 
locity, temperature, and density are quite consistent with 
the spectroscopic and interplanetary scintillation observa- 
tions (Fig. 2). Even though previous on e dimensional 
simulations (jSuzuki Sz Inutsukal I2005L 120061 ) show similar 
radial variations, the coronal heating and solar wind ac- 
celeration mechanism in our two dimensional simulation is 
essentially different from the previous ones. 

The energy losses such as radiative cooling, thermal 
conduction, and adiabatic cooling due to the solar wind 
are the main cooling processes in the solar atmosphere. 
In order to maintain the solar corona, heating processes 
are necessary to balance the cooling processes. As shown 
in the panel (a) of figure 3, the energy flux of the Alfven 
wave decreases monotonically, a part of which is trans- 
ferred to the solar wind. Although a sizable fraction of 
flux is decreased in the chromosphere by the reflection and 
dissipation processes, the energy flux that is supplied to 
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Fig. 2. — Comparison of the simulation and the ob- 
servation is summarized below. Red solid lines in each 
panel represent results of the simulations; these values 
are averaged over distance and over 30 minutes, (a) 
The gre en crosses dTeriaca est al.l I2003D and the blue tri- 
angles ([Zangrilli et al.l I2002D represent the proton out- 
flow speeds in the polar region observ ed by SOHO The 
black circles with crossed error bars (IGrall et al.l fl99l 
are obtai ned by VLBA. The b lack circles with vertical er- 
ror bars (|Habbal et al.l I1995D indicate m easurements by 
SPAR TAN 201-01. (b) The blue circles dBaneriee et all 
119981 ) show the nonthermal broa dening inferred f rom 
SUMER/SOHO. The black crosses (jEsser et al.lfl99l are 
derived empirically from nonthermal broadening based 
on t he UVCS/SOHO m easurements, (c) The black cir- 
cles (jFludra et al.l Il999f ) show elec tron temperature by 
CDS/SOHO. (d) Th e black circles dWilhelm et al.l ll998h 
and the blue stars dTeriaca et al.l I2003D are data based 
on observations by SUMER/SOHO and by CDS/SOHO 
respectively. The g reen triangles dTeriaca et all 120031) 
and orange squares (|Lamv et al.l 119971) are observed by 
LASCO/SOHO. 
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Fig. 3. — Alfven wave energy flux and its dissipation 
processes, (a) Alfven wave energy flux as a function of ra- 
dius, (b) Heating and cooling rate as a function of radius. 
The green solid line shows shock-heating rate estimated 
by counting sudden entropy jumps. The red-hatched area 
indicates turbulent heating rate estimated from Fourier 
and dimensional analysis. The black solid line and the 
blue solid line show radiative and adiabatic cooling rates, 
respectively. The black-hatched area represents the tran- 
sition region. 

the corona is well above 10 5 erg cm -2 s _1 , a typical num- 
ber t hat is required to maint ain the corona and the solar 
wind (jHansteen fc Leerlll995| ) . The dissipation mechanism 
of the Alfven wave should be different in various regions 
depending on the background medium. 

The panel (b) of figure 3 shows comparison of each com- 
ponent of the heating and cooling rates. The heating is 
separated into the compressive (dilatational) and incom- 
pressible (shearing) parts. The green solid line shows the 
compressive component; compressive waves are generated 
due to the nonlinear mode conversion fro m the Alfven wave 
to co mpressive (or slow mode) wave (|Kudoh fc Shibatal 
Il999t ) and these compressive waves eventually steepen into 
shocks. We estimate the heating rate by counting the 
entropy jumps at shock fronts in the simulations. The 
red-hatched area shows the incompressible heating, which 
is done by the dissipation of Alfvenic turbulence owing 
to strong shearing motion at small scales. We estimate 
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Fig. 4. — Alfven wave nonlinearity as a function of radius. 
Alfv'en wave nonlinearity is determined as (dva+)=(Ve — 
Be I \JATip)/2 1 and Va is Alfven speed. The black-hatched 
area represents the transition region. 

the turbulent heating rate from the Fourier component of 
Alfvenic disturbances. 

In the chromosphere, both shock and turbulence con- 
tribute to the heating. The wave nonlinearity, which is de- 
fined as wave amplitude divided by phase speed, quickly 
increases in the chromosphere with the rapid expansion 
of the magnetic flux tube (Fig. 4). As a result, com- 
pressive waves are generated effectively by the nonlinear 
mode conversion of Alfven waves. The turbulent heating is 
also important in the chromosphere because the Alfvenic 
turb ulence is developed efficien tly due to both phase mix- 
i ng dHevvaerts fc Priest! [T983|) as well as wave reflection 
(jMatthaeus et al.lll999t ). Since the flux tube rapidly opens 
near the photosphere, the Alfven speeds of the neighbor- 
ing field lines arc different with each other. Due to the 
difference in Alfven speeds across the magnetic field, the 
Alfven waves along the neighboring field lines gradually 
become out of phase, even though the waves are in phase 
at the photosphere, which creates strong shear to dissi- 
pate their wave energy. In addition to the phase mixing, 
the rapid decrease of the density in the chromosphere and 
the transition region causes increase in the Alfven speed 
that finally leads to the reflection of the Alfven wave. The 
nonlinear wave-wave interaction between the pre-existing 
outward component and the reflected component develops 
turbulent cascade. 

In the corona, the shock dissipation is the main con- 
tributor to the heating although the turbulence is also 
effective in the lower corona. Passing through the tran- 
sition region, the wave nonlinearity decreases rapidly be- 
cause of the large Alfven speed in the corona, so the local 
shock formation in the corona is not significant. Instead, 
compressive waves are generated by the vertically fluctu- 
ating motion of the transition region; the nonlinearly ex- 
cites longitudinal waves in the chromosphere conti nuously 
tap the transition region dKudoh fc Shibatal Il999h . which 
further excites upward propagating compressive waves in 
the corona. The reflected wave component drops off sig- 
nificantly above the transition region (Fig. 5) because 
the Alfven speed does not change so much. As a result, 
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Fig. 5. — Energy ratio of downward propagating Alfvcn 
wave to upward propagating Alfvcn wave. The energy of 
downward(+)/upward(-) propagating Alfvcn wave is pro- 
portional to (dva±) 2 = (Vfl =F B e / 'V&rp) / '4. The black- 
hatched area represents the transition region. 

the turbulent cascade is suppressed in the subsonic region 
(1.02 to 4 Rs). The power index of the energy spectral 
density of Alfvenic disturbance is significantly softer than 
-5/3, which indicates that the energy cascading to smaller 
scales is not effective in this region (Fig. 6). The shock 
heating compensates for the absence of turbulent heating 
in order to balance the cooling there. Generally, heating 
below the sonic point controls mass loading to the solar 
wind, so we suggest that the shock heating mechanism 
works efficiently to determine the mass loss rate from the 
sun. In our simulation, mass loss rate is of the order of 
10 -14 M0 yr _1 which agrees reasonably well with the ob- 
served value. 

In the solar wind acceleration region (R > 4 R s ) , again, 
both turbulent heating and shock heating are compara- 
ble. The wave nonlinearity, once dropped above the tran- 
sition region, increases gradually with radius due to the 
decrease in Alfven speed. Even though the nonlinearity 
is still small, Alfven waves suffer nonlinear effects due to 
their long traveling distance. Then, compressive waves are 
locally excited by the mode conversion, and these waves 
finally get dissipated by the shocks. The turbulent heat- 
ing works effectively in the solar wind. The signature of 
turbulent cascade can be seen clearly in the power spec- 
tra of Alfven waves (Fig. 7). Moreover, the turbulent 
cascade is not isotropic but an isotropic (jShebalin et al.l 
119831: iGoldreich fe Sridharl Il995l ) : the direction of turbu- 
lent cascading is perpendicular to the background mag- 
netic held. The energy cascading is triggered by the non- 
linear wave- wa ve interaction between the outgoing and re- 
flection waves (jMatthaeus et al.1ll999f ) . The increase of the 
reflection component of Alfven waves can be seen in figure 
7. The work done by the gas and the magnetic pressure 
from the turbulence is also of the same order and found 
to be sufficient to accelerate the solar wind. Figures 5 and 
7 indicate that the crit ical balance state, k r v a ~ k^vg, 
of Alfvenic turbulence (jGoldreich fe Sridharl Il995[ ) is at- 
tained as R increases from the wave-like state, or weak 
turbulence state, (k r VA > k^vg, at most) in the low corona 
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Fig. 6. — (a) Normalized power spectral density with 
respect to normalized the horizontal wave number. The 
symbols of plus, asterisk, diamond, and triangle show the 
power spectral density at 0.01, 0.1, 1, and 10 Rs, respec- 
tively, (b) Power index of the turbulent power spectral 
density with respect to height. The vertical axis represents 
the power index of Alfvenic turbulent spectrum. The hori- 
zontal dotted line indicates the values of Kolmogorov type 
turbulence, —5/3. The black-hatched area represents the 
transition region. 

fe.g.. lPerez fe Boldvrevl[2008fl . 

Our two-dimensional MHD simulation shows that both 
the shock and the turbulence are important for the coro- 
nal heating and the solar wind acceleration. We showed 
that the energy exchange between the Alfven mode and 
slow mode is effective, although previous MHD simulations 
of turbulence in homogeneous me dia show the decouplin g 
of Alfven and compressive modes (|Cho fe Lazaria n 2003). 
The inhomogeneity of the background medium due to the 
density stratification and the rapidly expanding flux tube 
is essential to understand the energy transfer processes 
in the solar atmosphere. In the previ ous ID MHD simula- 
tions (jSuzuki fc Inutsukall2005l |2006( ) , the shock heating is 
over-estimated because the geometrical expansion dilutes 
the shocks in a multidimensional system. We showed that 
the shock dilution is not so significant in a 2D system. 
The turbulent cascading process in our simulation results 
from 2D nonlinear terms while the previous studies (e.g. 
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Fig. 7. — Evolution of power spectrum of Alfven wave as 
a function of radius. The four panels represent the power 
spectrum of Alfven wave at (a) R=l.l Rs, (b) R=6.0 Rs, 
(c) R=ll Rs, and (d) R=14 Rs. The vertical and the 
horizontal axis represents the wave number in (j> direction, 
kj, , and in r direction, k r , normalized by fc = 2nRA(p, 
where Acf> is the angular system size in <f> direction. The 
color in each panel shows the power spectral density in 
Fourier space normalized by the peak value. 



iDmitruk fe Matthaeusl[2003HVan Ballegooijen et al.ll2011ft 

show the importance of nonlinear terms originating in 3D 
nature. Therefore 3D MHD simulation is necessary not 
only to completely rule out the possibility of shock dilu- 
tion but also to verify the turbulent cascading processes. 

We showed that these changes in the main contribu- 
tor to the heating are the natural consequences of the 
strong inhomogeneity, originated both in the photospheric 
magnetic field and in the gravitational stratification. Al- 
though only the Alfvenic disturbances are excited at the 
photosphere, in our simulation, compressive waves are 
also important especially for the chromospheric heatin g 
(|Bogdan et alJl2003t lHasan et al.ll2005b iFedun et al.ll2011l ) . 
It is possible to change the reflection rate of Alfven wave 
because the extra heating by compressive waves modifies 
the density structure around the transition region. 

Numerical computations were carried out on Cray XT4 
at Center for Computational Astrophysics, CfCA, of Na- 
tional Astronomical Observatory of Japan. Takuma Mat- 
sumoto gratefully acknowledges the research support in 
the form of fellowship from the Japan Society for the Pro- 
motion of Science for Young Scientists. 
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